{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Visualization of Coincident ICESat-2, Sentinel, and Argo Float Data\n",
    "\n",
    "ICESat-2 hackweek final project \n",
    "June 16, 2020  \n",
    "Nicole Abib"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "from icepyx import icesat2data as ipd\n",
    "import rasterio\n",
    "from rasterio.plot import show\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import h5py \n",
    "import sys  \n",
    "sys.path.insert(0, '/home/jovyan/leading_to_phytoplankton/scripts')\n",
    "import readers as rd\n",
    "import cartopy.crs as ccrs\n",
    "from astropy.time import Time\n",
    "from osgeo import gdal\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "# magic function to enable interactive plotting\n",
    "%matplotlib widget "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define path to sample data\n",
    "data_loc='/home/jovyan/shared/leading_to_phytoplankton/'\n",
    "argo_arc_fn = data_loc+'argo_arc.csv'\n",
    "argo_ant_fn = data_loc+'argo_ant.csv'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "NaN,NaN,2018,6,20,17,45,4.33,0.00079445,28.762,34.726,0.2117\n",
      "NaN,NaN,2018,6,20,17,45,5.6,0.00081808,28.764,34.726,0.2117\n",
      "NaN,NaN,2018,6,20,17,45,11.36,0.00081808,28.771,34.726,0.2117\n",
      "NaN,NaN,2018,6,20,17,45,16.69,0.00095985,28.773,34.726,0.2263\n",
      "NaN,NaN,2018,6,20,17,45,21.73,0.00084171,28.774,34.727,0.2044\n",
      "NaN,NaN,2018,6,20,17,45,26.77,0.00086534,28.774,34.727,0.219\n",
      "NaN,NaN,2018,6,20,17,45,31.73,0.00085352,28.744,34.728,0.2555\n",
      "NaN,NaN,2018,6,20,17,45,36.28,0.00087715,28.715,34.729,0.3139\n",
      "NaN,NaN,2018,6,20,17,45,41.36,0.00086534,28.702,34.727,0.3504\n",
      "NaN,NaN,2018,6,20,17,45,46.71,0.00084171,28.698,34.731,0.3942\n"
     ]
    }
   ],
   "source": [
    "# Quick check of csv file contents\n",
    "!head $argo_arc_fn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-50.01,131.59,2018,9,14,11,17,0.44,-0.00043419,7.505,34.27,-0.00013273\n",
      "-50.01,131.59,2018,9,14,11,17,1.54,-0.00042907,7.5033,34.267,-0.0035173\n",
      "-50.01,131.59,2018,9,14,11,17,2.6,-0.0004344,7.4997,34.271,0\n",
      "-50.01,131.59,2018,9,14,11,17,3.6,-0.00042339,7.493,34.278,0\n",
      "-50.01,131.59,2018,9,14,11,17,4.5,-0.00042338,7.495,34.276,-0.0073\n",
      "-50.01,131.59,2018,9,14,11,17,5.5,-0.00042338,7.4962,34.275,0\n",
      "-50.01,131.59,2018,9,14,11,17,6.5,-0.00042338,7.497,34.276,0\n",
      "-50.01,131.59,2018,9,14,11,17,7.4,-0.00043442,7.498,34.278,0\n",
      "-50.01,131.59,2018,9,14,11,17,8.54,-0.00043393,7.4952,34.277,-0.0069756\n",
      "-50.01,131.59,2018,9,14,11,17,9.54,-0.00042371,7.5001,34.277,-0.000146\n"
     ]
    }
   ],
   "source": [
    "!head $argo_ant_fn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "   lat  lon  year  month  day  hour  minute  Depth (m)  bbp (700 nm)  \\\n",
      "0  NaN  NaN  2018      6   20    17      45       5.60      0.000818   \n",
      "1  NaN  NaN  2018      6   20    17      45      11.36      0.000818   \n",
      "2  NaN  NaN  2018      6   20    17      45      16.69      0.000960   \n",
      "3  NaN  NaN  2018      6   20    17      45      21.73      0.000842   \n",
      "4  NaN  NaN  2018      6   20    17      45      26.77      0.000865   \n",
      "\n",
      "   temperature  salinity  chlorophyll  \n",
      "0       28.764    34.726       0.2117  \n",
      "1       28.771    34.726       0.2117  \n",
      "2       28.773    34.726       0.2263  \n",
      "3       28.774    34.727       0.2044  \n",
      "4       28.774    34.727       0.2190  \n",
      "     lat     lon  year  month  day  hour  minute  Depth (m)  bbp (700 nm)  \\\n",
      "0 -50.01  131.59  2018      9   14    11      17       1.54     -0.000429   \n",
      "1 -50.01  131.59  2018      9   14    11      17       2.60     -0.000434   \n",
      "2 -50.01  131.59  2018      9   14    11      17       3.60     -0.000423   \n",
      "3 -50.01  131.59  2018      9   14    11      17       4.50     -0.000423   \n",
      "4 -50.01  131.59  2018      9   14    11      17       5.50     -0.000423   \n",
      "\n",
      "   temperature  salinity  chlorophyll  \n",
      "0       7.5033    34.267    -0.003517  \n",
      "1       7.4997    34.271     0.000000  \n",
      "2       7.4930    34.278     0.000000  \n",
      "3       7.4950    34.276    -0.007300  \n",
      "4       7.4962    34.275     0.000000  \n"
     ]
    }
   ],
   "source": [
    "# Load the csv file with Pandas\n",
    "argo_arc_df = pd.read_csv(argo_arc_fn)\n",
    "argo_ant_df = pd.read_csv(argo_ant_fn)\n",
    "# Add column names defined in the metadata\n",
    "argo_arc_df.columns = ['lat', 'lon', 'year', 'month', 'day', 'hour', \n",
    "                   'minute', 'Depth (m)', 'bbp (700 nm)', 'temperature', \n",
    "                   'salinity', 'chlorophyll']\n",
    "argo_ant_df.columns = ['lat', 'lon', 'year', 'month', 'day', 'hour', \n",
    "                   'minute', 'Depth (m)', 'bbp (700 nm)', 'temperature', \n",
    "                   'salinity', 'chlorophyll']\n",
    "# Inspect the dataframe\n",
    "print(argo_arc_df.head())\n",
    "print(argo_ant_df.head())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# argo_df_nan = argo_df.replace(99999, np.nan)\n",
    "# argo_df = argo_df_nan"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "3f91914f63ee4aa88ab370fd35790f16",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Create a scatter plot showing data locations\n",
    "plt.figure(figsize=(8,8), dpi= 90)\n",
    "ax = plt.axes(projection=ccrs.PlateCarree()) # choose polar sterographic for projection\n",
    "ax.coastlines(resolution='50m', color='black', linewidth=1)\n",
    "plt.scatter(argo_arc_df['lon'], argo_arc_df['lat'], c= 'r',s=1,transform=ccrs.PlateCarree())\n",
    "plt.scatter(argo_ant_df['lon'], argo_ant_df['lat'], c= 'r',s=1,transform=ccrs.PlateCarree())\n",
    "ax.set_title('Location of Argo Floats');\n",
    "# plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Load in ICESat-2 Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n"
     ]
    }
   ],
   "source": [
    "is2_fn = 'IS2_S2/ATL03_20190805215948_05930404_002_02.h5'\n",
    "is2_f = h5py.File(data_loc+is2_fn, 'r')\n",
    "print(is2_f['orbit_info/sc_orient'][0]) # Check forward/backward orientation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Read in icesat 2 data as a pandas dataframe\n",
    "beam='gt1l'\n",
    "df03_b1=rd.getATL03(is2_f,beam)\n",
    "df03_b1.head()\n",
    "beam='gt2l'\n",
    "df03_b2=rd.getATL03(is2_f,beam)\n",
    "beam='gt3l'\n",
    "df03_b3=rd.getATL03(is2_f,beam)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load in Sentinel 2 Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "s2_fn = 'IS2_S2/2019085_T11XMK.tif'\n",
    "s2_img = rasterio.open(data_loc+s2_fn)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "# show(s2_img)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BoundingBox(left=399960.0, bottom=8890200.0, right=509760.0, top=9000000.0)\n",
      "EPSG:32611\n"
     ]
    }
   ],
   "source": [
    "s2_bounds  = s2_img.bounds\n",
    "print(s2_bounds)\n",
    "print(s2_img.crs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "POLYGON ((509760 8890200, 509760 9000000, 399960 9000000, 399960 8890200, 509760 8890200))\n"
     ]
    }
   ],
   "source": [
    "from shapely.geometry import box\n",
    "geom = box(*s2_bounds)\n",
    "print(geom.wkt)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "9521398d5f8d4f7f895b920333fe3503",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "var= 'heights' #choose which variable we want to plot\n",
    "\n",
    "## we will want to set colorbar parameters based on the chosen variable\n",
    "vmin=-10\n",
    "vmax=30\n",
    "ticks=np.arange(-20,100,5)\n",
    "\n",
    "plt.figure(figsize=(8,8), dpi= 90)\n",
    "ax = plt.axes(projection=ccrs.NorthPolarStereo(central_longitude=0)) # choose polar sterographic for projection\n",
    "ax.coastlines(resolution='50m', color='black', linewidth=1)\n",
    "ax.set_extent([-180, 180, 60, 90], ccrs.PlateCarree())\n",
    "plt.scatter(argo_arc_df['lon'], argo_arc_df['lat'], c= 'r',s=2,transform=ccrs.PlateCarree())\n",
    "plt.scatter(df03_b1['lons'][::1000], df03_b1['lats'][::1000],c=df03_b1[var][::1000], cmap='viridis', vmin=vmin,vmax=vmax,transform=ccrs.PlateCarree())\n",
    "plt.plot(*geom.exterior.xy, color='black', linewidth=1,transform=ccrs.UTM(11))\n",
    "plt.colorbar(label=var, shrink=0.5, ticks=ticks,extend='both');"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "ac2a551da25c4249893e2dffbf24bf18",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "var= 'heights' #choose which variable we want to plot\n",
    "\n",
    "## we will want to set colorbar parameters based on the chosen variable\n",
    "vmin=-10\n",
    "vmax=30\n",
    "ticks=np.arange(-20,100,5)\n",
    "\n",
    "# plt.figure(figsize=(8,8), dpi= 90)\n",
    "fig,ax = plt.subplots(figsize=(8,8), dpi=90)\n",
    "ax = plt.axes(projection=ccrs.NorthPolarStereo(central_longitude=0)) # choose polar sterographic for projection\n",
    "\n",
    "ax.coastlines(resolution='50m', color='black', linewidth=1)\n",
    "# ax.set_extent([-180, 180, 60, 90], ccrs.PlateCarree())\n",
    "plt.scatter(argo_arc_df['lon'], argo_arc_df['lat'], c= 'r',s=2,transform=ccrs.PlateCarree())\n",
    "plt.scatter(df03_b1['lons'][::1000], df03_b1['lats'][::1000],c=df03_b1[var][::1000], cmap='viridis', s=2, vmin=vmin,vmax=vmax,transform=ccrs.PlateCarree())\n",
    "plt.scatter(df03_b2['lons'][::1000], df03_b2['lats'][::1000],c=df03_b2[var][::1000], cmap='viridis', s=2, vmin=vmin,vmax=vmax,transform=ccrs.PlateCarree())\n",
    "plt.scatter(df03_b3['lons'][::1000], df03_b3['lats'][::1000],c=df03_b3[var][::1000], cmap='viridis', s=2, vmin=vmin,vmax=vmax,transform=ccrs.PlateCarree())\n",
    "plt.plot(*geom.exterior.xy, color='black', linewidth=1,transform=ccrs.UTM(11))\n",
    "ax.set_extent([ -125, -115, 80, 81])\n",
    "plt.colorbar(label=var, shrink=0.5, ticks=ticks,extend='both');"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
